Introduction

Deforestation is accelerating at an alarming rate around the world as a result of Palm oil production. In order to satisfy growing global demand for Palm oil, an ingredient commonly found in cosmetics to processed foods - clearing the way for lucrative oil plantations by forest fires is negatively impacting our environment. Concerns include the loss of rainforest habitat for endangered species, harmful carbon dioxide emissions, coupled with the displacement of local indigenous communities.

Exploratory Questions

What countries are the top producers of Palm oil? Where are the palm oil mills located? Where are the forest-fire hotspots?

The Datasets

The UN Food and Agriculture organization provides information on global Palm oil production and area harvested. The Global Forest Watch also has Palm oil Mill location data and fire alert data available for our review.

Industry Overview

To navigate to where Palm Oil production is at its highest volume, we first approach the UN Food and Agriculture dataset:

setwd("/Users/drscholls303/Desktop/Project1")
# Load the required R packages and specified datasets
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(ggplot2)
## Warning: package 'ggplot2' was built under R version 3.2.4
faostat.palm <-read.csv("Faostats.csv")
faostat.palm <- faostat.palm %>% na.omit()
#Select relevant columns
producers <-select(faostat.palm, AreaName, ElementName, Value, Year)
# Filter for top producers in most recent available year of 2014
producers.year<-filter(producers, Year == 2014) %>% 
  group_by(AreaName, ElementName='Production') %>% 
  summarise(., Production=sum(Value)) 
nrow(producers.year) #46 rows 
## [1] 46
#Calculate percentage of world production for each country
producers.year$Percent = round(producers.year$Production[1:46]*100/
          sum(producers.year$Production), digits = 2)
producers.year<-tbl_df(producers.year)
top<-top_n(producers.year,5,Production)
top
## Source: local data frame [5 x 4]
## 
##    AreaName ElementName Production Percent
##      (fctr)       (chr)      (int)   (dbl)
## 1  Colombia  Production    6998364    1.99
## 2 Indonesia  Production  163514286   46.52
## 3  Malaysia  Production  120627960   34.32
## 4   Nigeria  Production   11020724    3.14
## 5  Thailand  Production   13355542    3.80
#top 5 graph
library(scales)
top.plot<-ggplot(top, aes(AreaName, Production, fill=Percent)) + geom_bar(stat="identity") + xlab("Country") + ylab("Oil Production") +
  ggtitle("Palm Oil Production by Country (Top 5)")
  
firstplot<- top.plot + 
scale_y_log10(
breaks = scales::trans_breaks("log10", function(x) 10^x),
labels = scales::trans_format("log10", math_format(10^.x)))
firstplot

The top five producers in the world include Indonesia and Malaysia. Indonesia accounts for over 46% or worldwide production. Malaysia accounts for around 34%. The combined 80% allocation in this region should be noted.

Industry Map - Bird’s eye view

library(maptools)
## Loading required package: sp
## Checking rgeos availability: TRUE
library(ggplot2)
library(ggmap)
library(rgeos)
## rgeos version: 0.3-11, (SVN revision 479)
##  GEOS runtime version: 3.4.2-CAPI-1.8.2 r3921 
##  Linking to sp version: 1.1-0 
##  Polygon checking: TRUE
library(dplyr)
plotmap <- readShapePoly(fn="TM_WORLD_BORDERS-0.3/TM_WORLD_BORDERS-0.3.shp")
palmHa <-read.csv("Faostats.csv", stringsAsFactors=F)
palmHa <-filter(palmHa, ElementName == 'Area harvested')
names(plotmap)[5] <- "AreaName"
target <- c("Indonesia", "Malaysia", "Colombia", "Thailand", 
            "Philippines", "Nigeria", "Cameroon", "Ecuador", "Ghana", "Honduras","Guatemala")
#Use foritfy to transform spatial plotmap into dataframe
plotmapDf <- fortify(plotmap, region = "AreaName") %>% 
  filter(id %in% target)
palmHaMapDf <- merge(plotmapDf, palmHa, by.x="id", by.y="AreaName")
backgroundMap <- ggplot(data=palmHaMapDf) + geom_path(aes(x=long, y=lat, group=group), color='grey') + coord_equal() + geom_polygon(aes(x=long, y=lat, group=group, fill=Value))
mapAsia <- get_map(location = 'Indonesia', zoom=4)
## Map from URL : http://maps.googleapis.com/maps/api/staticmap?center=Indonesia&zoom=4&size=640x640&scale=2&maptype=terrain&language=en-EN&sensor=false
## Information from URL : http://maps.googleapis.com/maps/api/geocode/json?address=Indonesia&sensor=false
mapAfrica <- get_map(location = 'Africa', zoom=4)
## Map from URL : http://maps.googleapis.com/maps/api/staticmap?center=Africa&zoom=4&size=640x640&scale=2&maptype=terrain&language=en-EN&sensor=false
## Information from URL : http://maps.googleapis.com/maps/api/geocode/json?address=Africa&sensor=false
mapSouthAm <- get_map(location = 'South America', zoom=4)
## Map from URL : http://maps.googleapis.com/maps/api/staticmap?center=South+America&zoom=4&size=640x640&scale=2&maptype=terrain&language=en-EN&sensor=false
## Information from URL : http://maps.googleapis.com/maps/api/geocode/json?address=South%20America&sensor=false
ggmapObjAsia <- ggmap(mapAsia)
Asia <- ggmapObjAsia + geom_polygon(aes(x=long, y=lat, group=group, fill=Value), data=palmHaMapDf, alpha=.9) + 
  geom_path(aes(x=long, y=lat, group=group), data=palmHaMapDf, color='black')
Asia

ggmapObjAf <- ggmap(mapAfrica)
Africa <- ggmapObjAf + geom_polygon(aes(x=long, y=lat, group=group, fill=Value), data=palmHaMapDf, alpha=.9) + 
  geom_path(aes(x=long, y=lat, group=group), data=palmHaMapDf, color='black')
Africa
## Warning: Removed 22194 rows containing missing values (geom_path).

ggmapSouthAm <- ggmap(mapSouthAm)
SouthAm <- ggmapSouthAm + geom_polygon(aes(x=long, y=lat, group=group, fill=Value), data=palmHaMapDf, alpha=.9) + 
  geom_path(aes(x=long, y=lat, group=group), data=palmHaMapDf, color='black')
SouthAm
## Warning: Removed 30888 rows containing missing values (geom_path).

As a visual confirmation on the map - looking at SouthEast Asia vs Africa or South America, there are more acres harvested within Indonesia.

Palm Oil Harvested in Indonesia

library(dplyr)
library(ggplot2)
faostat.palm <-read.csv("Faostats.csv")
faostat.palm <- faostat.palm %>% na.omit()
AH<-filter(faostat.palm, ElementName=='Area harvested', AreaName=='Indonesia')

areaGraph<-ggplot(AH, aes(Year, Value, fill=AreaName)) + geom_bar(stat="identity") +
  ylab("Area Harvested") + ggtitle("Area Harvested through Time(Indonesia)")
areaGraph

The rise in Area Harvested over time is alarming over time. The drastic percent change from 2000 to 2014 was over 186%. From a value of 2,014,000 from 740,709.

Mill Distribution and Fire Hotspots

library(dplyr)
library(leaflet)
mills<-read.csv("mills.csv")
fires<-read.csv("fires.csv")
last.year.fires <-filter(fires, OBJECTID > 11556)#fire-alerts, as of 1/01/2015
#zoom in on SE asia, Circles for mills, Markers for fire alerts
map <- leaflet(mills) %>%
  addTiles('http://{s}.basemaps.cartocdn.com/dark_all/{z}/{x}/{y}.png',attribution='Map tiles by <a href="http://stamen.com">Stamen Design</a>, <a href="http://creativecommons.org/licenses/by/3.0">CC BY 3.0</a> &mdash; Map data &copy; <a href="http://www.openstreetmap.org/copyright">OpenStreetMap</a>') %>% addCircles(~longitude, ~latitude, popup=mills$type, weight = 3, radius=40, color="#ffa500", stroke = TRUE, fillOpacity = 0.8) %>%
  addMarkers(data = last.year.fires, ~Latitude, ~Longitude, clusterOptions = markerClusterOptions()) %>% 
  addLegend("bottomright", colors= "#ffa500", labels="Mills", title="Legend") 
map

Conclusion

There is a strong presence of fire hotspots within the last year where the majority of Mills are located. Indonesia possesses the most Palm Oil Mills as currently reported in available data. There are hardly any visible fire hotspots in non-Palm Oil producing countries.

Future Questions

Are there any geographic overlaps with existing nature reservations- what is the tree cover loss over time in those areas? How much has the Air Quality Index been impacted? Is there a trend between species extinction and habitat location in Palm Oil producing nations?